function heat_step

%  plots solution of the heat equation, including surface and time slices
%       diff(u,x,x) = diff(u,t)   for xL < x < xR, 0 < t < tmax
%  where
%      u = 0  at x=xL,xR  and  u = step(a,b)  at t = 0

% set parameters
N=20;
M=20;
tmax=0.1;
xL=0;
xR=1;
a=0.25;
b=0.75;

% generate the points along the x-axis, x(1)=xL and x(N)=xR
x=linspace(xL,xR,N);
h=x(2)-x(1);

% generate the points along the t-axis, t(1)=0 and t(M)=tmax
t=linspace(0,tmax,M);
k=t(2)-t(1);

% evaluate the solution
U=zeros(N,M);
for i=1:N
    U(i,1)=0.5*(sign(x(i)-a)-sign(x(i)-b));
end

for j=2:M
    for i=1:N
        s=0;
        for ii=1:50
            an=2*(cos(a*pi*ii)-cos(b*pi*ii))/(pi*ii);
            s=s+an*exp(-pi*pi*ii*ii*t(j))*sin(pi*ii*x(i));
        end
        U(i,j)=s;
    end
end

% get(gcf)
set(gcf,'Position', [753 616 560 530]);
%subplot(2,1,1)
axes('position',[.13  .43  .8  .53])
surf(x,t,U')
%surfc(x,t,U')
%contour(x,t,U')
xlabel('x-axis','FontSize',14,'FontWeight','bold')
ylabel('t-axis','FontSize',14,'FontWeight','bold')
zlabel('Solution','FontSize',14,'FontWeight','bold')
% Set the fontsize to 12 for the plot '
set(gca,'FontSize',12);

% calculate the time slice solutions
% set parameters
N=100;

% pick time points (by picking the index)
t(1)=0.001;
t(2)=0.01;
t(3)=0.1;

% plot results
x=linspace(xL,xR,N);

%subplot(2,1,2)
axes('position',[.13  0.08  .8  .25])
hold on

xlabel('x-axis','FontSize',14,'FontWeight','bold')
ylabel('Solution','FontSize',14,'FontWeight','bold')

for i=1:N
    u(i)=0;
    for ii=1:50
        an=2*(cos(a*pi*ii)-cos(b*pi*ii))/(pi*ii);
        u(i)=u(i)+an*exp(-pi*pi*ii*ii*t(1))*sin(pi*ii*x(i));
    end
end

plot(x,u,'-k')

for i=1:N
    u(i)=0;
    for ii=1:50
        an=2*(cos(a*pi*ii)-cos(b*pi*ii))/(pi*ii);
        u(i)=u(i)+an*exp(-pi*pi*ii*ii*t(2))*sin(pi*ii*x(i));
    end
end

plot(x,u,'--k')
for i=1:N
    u(i)=0;
    for ii=1:50
        an=2*(cos(a*pi*ii)-cos(b*pi*ii))/(pi*ii);
        u(i)=u(i)+an*exp(-pi*pi*ii*ii*t(3))*sin(pi*ii*x(i));
    end
end

plot(x,u,'-.k')
box on
% Set the fontsize to 14 for the plot
set(gca,'FontSize',14);
legend(' t = 0.001',' t = 0.01',' t = 0.1',2);
% Set legend font to 14/bold
set(findobj(gcf,'tag','legend'),'FontSize',14,'FontWeight','bold');

hold off
